La ville de Paris rend publiquement accessible de vastes données sur la circulation routière dans la ville. La découverte du site web Paris Data (“Paris Data” visité le 26.12.2021) nous a motivé pour utiliser cette ressource dans notre projet de prédiction avec Yanig Goude. Ce qui nous a surtout intrigué est la possibilité de relier les différents capteurs dans une structure de graphe. D’un côté, cela motive des études globales du réseau entier : Est-ce qu’on trouve un algorithme global, qui apprend les interactions dans le réseau entier ? De l’autre côté, l’examination de rues individuelles mène à des questions sur la corrélation locale de la circulation. Un grand boulevard se comporte-t-il comme ses avenues voisines ? Enfin, dans tout cela, nous nous poserons la question à quel horizon temporel nous réussirons à prédire la présence de voitures dans les rues de Paris.
Avant de nous lancer, nous présentons d’abord le jeu de données en détail, sa richesse mais aussi les difficultés qu’il apporte : Sa taille et les valeurs manquantes. Ayant passé une grande partie du projet sur cet aspect, nous présentons nos stratégies pour résumer et compléter ce data set.
Cela fait, nous formulons plusieurs problèmes de prédiction, à différents horizons temporels et en utilisons différentes parties du data set. (TODO: spécifier à la fin pour avoir une petite table de matières ici)
Sur le site Paris Data, on trouve énormité de jeux de données issus des activités de la ville, notamment sur le de comptage routier (“Données Du Comptage Routier” visité le 26.12.2021). Ces données collectées grâce à des boucles électromagnétiques implantées dans la chaussée sur plus de 3000 tronçons de voies. Ces boucles et bornes de comptage mesurent principalement deux choses: Le nombre de voitures qui passent et l’occupation du tronçon de route. Disponibles de 2014 jusqu’aujourd’hui, ces deux chiffres sont mesurés au pas horaire, résultant en environ \[ 8 \text{ ans } \times 365 \text{ jours } \times 24 \text{ heures } \times 3000 \text{ capteurs } = 2.1\times10^8 \text{ lignes dans le data set.} \]
Chacune de ces lignes contient les informations suivantes qui serons plus tard nos variables explicatives \(X_i\). Les noms de ces variables seront toujours notés en gras.
Les données historiques de circulation sont enregistrés sur le site web sous le format .txt. Chaque fichier couvre toutes les bornes pour une semaine. Pour commencer, il faut donc tout télécharger et convertir de manière automatisée ces fichiers en dataframes. En ce faisans, nous effectuons également un triage par identifiant de borne iu_ac.
En vue de la quantité énorme de 200 millions d’observations, nous sommes obligés à réduire la taille du dataset. La première démarche consiste à retirer toutes les variables inintéressantes, laissant t_1h, nbCar, rateCar et iu_ac.
Ayant réduit la dimension de chaque point dans les données, le deuxième point qu’on peut attaquer est le nombre gigantesque de senseurs de comptage. Pour ce faire, nous utilisons un modèle simplifié des rues de Paris. Comme le dataset couvre la plupart des grandes axes routières, y inclus le périphérique, notre objectif était d’agréger les données au long de ces rues :
Pour identifier où passent le plus de voitures, nous calculons des moyennes sur nbCar pour chaque libellé, c’est à dire la moyenne à travers le temps et toutes les bornes associées à chaque libellé. En prenant les 200 libellés les plus utilisées en termes de nombre de voitures (en excluant le périphérique), nous obtenons le graphique 2.1 qui représente très bien les axes routières auxquelles on s’attendrait.
Ensuite, nous ignorons les noms des libellés et regardons exclusivement ce graphique. De manière arbitraire nous en déduisons le graphe dans figure 2.2 qui représente alors un schéma très simplifié de la circulation à Paris. Grâce à notre démarche d’utiliser les rues les plus fréquentés, cette abstraction nous permet d’avoir un modèle représentatif de la circulation parisienne.
Figure 2.1: Représentation en rouge des bornes d’observation associés aux 200 libellés les plus fréquentés (hors périphérique)
Figure 2.2: Graphe simplifié. Points noirs indiquent les intersections entre nos arêtes.
Chaque arête dans le graphe 2.2 regroupe plusieurs senseurs dont nous agrégeons les données : Après une collecte minutieuse et laborieuse des identifiants de toutes les bornes sur chacune des arêtes, nous associons pour chaque heure la moyenne de nbCar et rateCar des bornes à l’arête dont ils font partie. Cela a surtout trois avantages:
Ayant suffisamment réduit les données pour qu’ils rentrent dans le RAM d’un PC, il reste pourtant un problème: La présence de données manquantes.
En examinant des plots basiques des taux d’occupation, nous voyons souvent une image comme dans le graphique TODO. Malgré l’agrégation au long des arêtes, il reste un nombre considérable de NA dans le dataset. Afin d’executer des tâches de prédiction, il nous faut pourtant des données complètes et une partie importante de notre travail était de remplir les trous.
Figure 2.3: TODO remplacer graphique par un plot du taux de circulation d’un jour avec un trou.
Figure 2.4: TODO remplacer graphique par distribution des valeurs manquantes pour voir les deux types de trous (todo : Jakob)
En regardant le graphique (fig:NAdistribution) qu’il y a deux types de “trous” dans le data set:
Le fait que les arêtes ont de différentes structures de valeurs manquantes peut aussi être observé en regardant la répartition du pourcentage de NA’s à travers les dataframes dans figure (fig:histograms).
Figure 2.5: Distribution des valeurs manquantes
Nous observons que les données de taux d’occupation sont généralement plus complètes : Ceci est dû au fait qu’il y a plus de capteurs de ce type. Par ailleurs, nous avons vérifié qu’une valeur manquante de nbCar implique toujours un NA dans rateCar.
Pour remplir les trous, la première étape consiste à rajouter des variables explicatives supplémentaires. En fonction du type de trou, nous verrons qu’ils auront un impact important sur la complétion.
A l’instant, nous disposons des variables nbCar et rateCar pour chaque arête et chaque heure entre 2014 et 2020. Nous exploitons deux propriétés structurelles des données pour obtenir de l’information supplémentaire:
D’abord, il s’agit de séries temporelles, donc nous rajoutons valeurs historiques de nbCar et rateCar décalés d’une heure, d’un jour et d’une semaine (nbCarLaggedHour, nbCarLaggedDay, et nbCarLaggedWeek) comme variables explicatives. L’espoir étant que la circulation se ressemble à travers des semaines consécutives par exemple. Avec le pure objectif de compléter les données manquantes, il est aussi pertinent de rajouter des valeurs du futur (rateCarFuturHour etc.), vu que le modèle d’imputation a pour but d’interpoler au lieu d’extrapoler. On a l’hypothèse que toutes ces variables décalées serviront surtout à remplir des trous locaux : Une heure ou une journée qui manque devrait facilement se laisser inférer en utilisant des données qui sont proches temporairement.
Pour les grands trous, nous n’avons souvent pas accès aux heures et jours d’avant. Cependant, grâce à l’interconnexion de nos données à travers le graphe présenté en haut, les mesures de circulation d’une arête peuvent bien être expliquées en fonction de celle de ces voisines. Dans un effort manuel, nous avons pour chaque arête \(A\) collectionné des voisins \(V_i\), dans le sens que les voitures dans \(A\) passent ensuite à l’un des \(V_i\) ou bien l’inverse. Nous faisons ici encore des choix arbitraires qui est voisin de qui, pour deux raisons. Premièrement, inférer statistiquement quelle arête influence quelle autre nécessiterait des données déjà complètes pour une régression par exemple. Deuxièmement, utiliser toutes les autres 68 arêtes au lieu de juste 2-6 voisins augmenterait trop le temps d’exécution de la complétion.
Donnons un exemple illustrant le choix de voisins ainsi que la dénomination des nouvelles variables : Parmi les voitures qui entrent dans Paris sur l’arête “pont amont - pont austerlitz” il y a une partie considérable qui continuent tout droit sur l’arête “pont austerlitz - chatelet.” Par conséquent, nous rajoutons la variable rateCar_pontausterlitzTOchatelet au dataframe de “pont amont - pont austerlitz.”
Jusqu’ici, nous n’avons utilisé que des informations déjà présentes dans le jeu de données. Dans la suite, nous rajoutons des variables extérieures qui pourraient également expliquer la circulation routière.
Variables temporelles
Le traffic routier étant relié à l’activité humaine, nous avons ajouté de nombreuses variables temporelles (notamment grâce à la librairie lubridate).
Index de la situation sanitaire en rapport avec le Covid-19
Nous avons également ajouté, covidIndex, un index allant de 0 à 100 représentant les restrictions du gouvernement sur la situation sanitaire en rapport avec le Covid-19. Il est calculé à partir de nombreux indicateurs et est fourni par l’université d’Oxford (“Covid Index” visité le 26.12.2021). Cette variable permettra éventuellement une étude de l’année 2020 qui voit une circulation perturbé. Pour cela, il faudra utiliser des modèles qui s’adaptent très vite à l’influence de nouvelles variables.
Météo
A partir de données de l’Organisation Météorologique Mondiale (“Observations Météorologiques Historiques En France” visité le 26.12.2021), nous avons extrait 2 variables météorologiques : temperature la température en Kelvin et precipitation les précipitations dans les 3 dernières heures en mm. Ces données ont été mesurées à Athis-Mons en Essonne.
On dispose d’un relevé tous les 3 heures environ donc on procède à une interpolation pour compléter les données. Etant donné que les mesures sont uniformément réparties, on utilise une interpolation linéaire basique à l’aide de la fonction na.approx de la librairie zoo.
TODO Jakob à rédiger avec graphiques (train dimanche).
Choix des jeux de données
On découpe notre jeu de données en 2 parties de proportion 2/3 et 1/3 : une partie apprentissage de 2014 à 2017 et une partie test de 2018 à 2019. Ainsi, on pourra évaluer les performances de nos modèles sur les différentes périodes d’une année.
Choix de la métrique
Pour évaluer ces performanes, on choisit la métrique Root Mean Square Error (RMSE) qui représente la racine carrée du deuxième moment d’échantillonnage des différences entre les valeurs prédites et les valeurs observées. On n’utilise pas la métrique Mean Absolute Percentage Error (MAPE) car un nombre important des valeurs de rateCar sont nulles.
Modèles naïfs à battre
Afin de comparer nos modèles, on construit 3 modèles témoins. Ces modèles sont dit naïfs car extrêmement simpliste :
naiveModel1 prévoit le taux d’occupation d’une heure \(t \in \{ 0, \dots, 23 \}\) d’un jour \(j \in \{ 1, \dots, 7\}\) d’un mois \(m \in \{ 1, \dots, 12 \}\) d’une année \(a \in \{ 2018, 2019 \}\) en moyennant les taux d’occupation de l’heure \(t\) du jour \(j\) du mois \(m\) pour toutes les années \(a \in \{ 2014, \dots, 2017 \}\).
naiveModel2 fait de même en calculant la médiane et non la moyenne.
naiveModel3 est simplement le taux d’occupation à l’heure précédente rateCarLaggedDay.
On calcule la moyenne des erreurs RMSE pour les 69 arêtes :
| naiveModel1 | naiveModel2 | naiveModel3 |
|---|---|---|
| 5.844 | 5.899 | 3.955 |
Recherche de paramètres
Afin d’optimiser les performances de nos différents modèles, nous effectuons de la recherche de paramètres à l’aide du package caret. On utilise les 2 méthodes suivantes :
cv : validation croisée sur 16 blocs. On divise ainsi nos 4 années de données d’apprentissage en blocs de 3 mois, ce qui englobe les différentes saisonnalités de nos données.
timeSlice : où on fixe une fenêtre initiale d’apprentissage de 2 ans, que l’on incrémente dans l’ordre chronologique de 2 mois à chaque itération, pour entraîner le modèle et mesurer sa performance à l’horizon 1 (c’est-à-dire prédire le taux d’occupation de l’heure suivante).
L’avantage de timeSlice est qu’il permet d’optimiser un modèle pour la prédiction à un horizon déterminé. De plus, contrairement à la validation croisée, il n’utilise pas de données futures pour prédire le passé. Cependant, en pratique, on souhaite entraîner le modèle une unique fois sans le mettre à jour tous les 2 mois. Dans ce cas, la validation croisée est plus pertinente.
Lors de nos recherches de paramètres, nous avons observé que les 2 méthodes sélectionnent des paramètres proches. Finalement, on se base sur la validation croisée en 16 blocs pour la recherche de paramètres.
Pour permettre la reproductibilité des résultats, on fixe la graine du générateur aléatoire. De plus, on effectue la recherche de paramètre pour une seule arête sur les 69 car c’est un processus long et on fait l’hypothèse que le paramètre optimal pour la prévision sur une arête est proche de celui de toutes les autres arêtes.
Figure 3.1: Recherche du paramètre cp optimal